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An effective Hartree-Fock-Bogoliubov-type interaction is calculated for the 
Hubbard model in second order in the coupling by means of flow equations. 
A stability analysis is performed in order to obtain the transition into various 
possible phases. 

We find, that the second order contribution weakens the tendency for the 
antiferromagnetic transition. Apart from a possible antiferromagnetic tran- 
sition the d-wave Pomeranchuk instability recently reported by Halboth and 
Metzner is usually the strongest instability. A newly found instability is of 
p-wave character and yields band- splitting. In the BCS-channel one obtains 
the strongest contribution for dx2_y2-waves. Other types of instabilities of 
comparable strength are also reported. 

PACS numbers: 71.10.Fd, 71.27. +a, 74.20-z. 



1. INTRODUCTION 

The Hubbard model is commonly considered a model which contains es- 
sential features of the Cuprate layers of high-temperature superconductorsJH 

In this paper we analyze the Hubbard model by means of flow equa- 
tion^. The basic idea is to eliminate off-diagonal matrix-elements of the 
interaction in order to obtain a diagonal or a block-diagonal Hamiltonian. 
This is done by a continuous unitary transformation of the Hamiltonian as 
a function of the flow-parameter I, 

^ = W),H{1)]. (1) 



In the original papem the Hamiltonian H was decomposed into two 
parts, the diagonal part and the off-diagonal part and the generator 
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rj of the flow equation was chosen so that H^- was ehminated. This could be 
done for example by choosing 

V=m,H]. (2) 

In the case of interacting electrons the diagonal terms in (P) were those 
conserving the number of quasi-particles, i.e. electrons above and holes 
below the Fermi-sea. Thus the contributions creating and annihilating quasi- 
particles across the Fermi-surface were eliminated. This scheme was useful 
at zero temperature. However, at finite temperatures the Fermi-surface is 
no longer sharp. Then this scheme is no longer useful. 

Instead we will choose a different approach. We rewrite eq.(§) 

V=[H,H,], (3) 

and in contrast to our previous approach we consider no longer terms to 
be either diagonal or off-diagonal, but we introduce a continuous quantity 
r > associated to the terms in the Hamiltonian Hj- which determines how 
urgently we wish to eliminate the corresponding term, r = means, that 
we keep the term. The larger r the more urgently we eliminate it. We call 
r the elimination factor. In our previous scheme the terms in the diagonal 
part -ffd have r = 0, the terms in the off-diagonal part Hj- have r = 1. From 
now on r may be continuous. 

We will use this scheme in order to transform the Hamiltonian into 
the form of a molecular-field type Hamiltonian, that is into a form which 
depends only on biquadratic terms c^Cfc, c|.cl^, and c_fcCfc, which yields a 
Hartree-Fock-Bogoliubov-type interaction. In order to accomodate antifer- 
romagnetism, too, we will keep terms c|._,_g^Cfc, and also c|.cj^_^, and Cgg-kCk, 
where qo = is the staggered wave-vector for lattice spacing a. This 

yields an effective interaction which contains the interaction terms described 
in eqs. (75,81,8^, and 95). In this way we eliminate the fiuctuations around 



the molecular-field type behavior. The calculation of this effective interac- 
tion is performed in second-order in the Hubbard-coupling U. 

This scheme differs from the renormalization-group scheme described by 
Shankari and applied by Zanchi and Schulzi, by Salmhofeil and by Halboth 
and MetznerH, and by Honerkamp et al.Q to the Hubbard model in so far, 
as they eliminate the degrees of freedom outside a shell enclosing the Fermi- 
surface. In our scheme however, first, the interactions connecting states far 
apart in energy, irrespective of their distance from the Fermi-level, are elimi- 
nated. This is in the same spirit as the similarity renormalization by Glazek 
and Wilson.l3 Since we transform the interaction we always deal with an effec- 
tive interaction and not with truncated and partially integrated correlations. 
This implies, that the interaction remains finite at phase-transitions. 
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Starting from the Hartree-Fock-Bogoliubov interaction we perform a 
stability analysis, that is we investigate, at which critical (U/t)c becomes the 
system unstable against any symmetry breaking contributions. It has long 
been known, that at half-filling the system becomes antiferromagnetic, and 
during recent years many calculations have shown, that the Hubbard model 
shows an attractive interaction for d^2_y2-wa,ve pairing. Apart from the 
above meri tiop edBJl We refer to the FLEX-approach by Bickers, Scalapino, 
and WhitaBll2l, the review by Scalapinoa and calculations by Hanke et alx3. 
We find that often the strongest instability (apart from antiferromagnetism) 
is a d-wave Pomeranchuk-instability recently reported by Halboth and Metz- 
neilll and a band-splittinginstability of p-wave character. The d-wave flux- 
phase instability (KotliarEJ, Affleck and MarstonEJ, Chakravarty et al.£3) 
appears degenerate with the superconductiving instability at half filling, and 
above at doping. 

In section 2 we derive a general expression for the generator r] from 
a gradient procedure and introduce the elimination factor r. In section 3 
the second-order contribution to the effective electron-electron interaction 
is derived. In section 4 the expressions for the free energy of the Hubbard 
model are given. In section 5 the numerical results for the phase-instabilities 
are presented. 



2. GENERATOR r] FOR FLOW EQUATIONS FROM A 
GRADIENT PROCEDURE 

In order to determine the generator rj in the flow equation (|^) we intro- 
duce a quadratic form of the Hamiltonian 

G{H) = ^Y.9i3,kiHjiHik. (4) 

This form is chosen in such a way that G becomes a minimum by means of 
the flow equations. Without restriction of the general procedure we require 
the symmetry 

9ij,kl = 9kl,ij- (5) 
G should be real. Since H is hermitean we have 

G*{H) = ^Y.9:,,klH^JHM. (6) 
Thus G* = G is obtained for 



9ij,kl — 9ji^lk- 



(7) 
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We consider now the variation of G under the flow 



dG 
"dT 




(8) 



Since dG/dl should be (semi-)negative, we choose 




which can be rewritten 




with 



Hr,ji — gij,lkHkl- 



(11) 



G should decrease under the flow. Therefore one chooses a quadratic form, 
which vanishes for the desired diagonal form, but whose other eigenvalues 
are positive. For this purpose one may use a double commutator with a 
hermitean operator v, 



With this choice those contributions to H are considered diagonal, which 
commute with v. If for example v is the operator of the quasi-numbers, then 
a matrix element Hij which connects states with vi and Vj quasi-particles 
contributes {vi — Vj)'^HijHji to G{H). For a general basis one has 



This matrix g is symmetric and hermitean as required. 

Instead of using one double-commutator one may use a sum of such 
contributions. This can be used for the investigation of systems at finite 
temperatures where due to the disappearance of a sharp Fermi-surface the 
concept of particles above and holes below the Fcrmi-surface is no longer 
valid. We describe the fermion system in terms of creation and annihilation 
operators and c and introduce the operator v as a one particle operator 



H^ = [v,[v,H]\. 



(12) 



9ji,ki = {v^)jk^ii - "^Vjuvii + Sjk{v^)ii. 



(13) 




(14) 




(15) 
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It indicates how strongly we wish to eliminate this term. The larger it is the 
more we consider it to be non-diagonal. If we now choose several of these 
contributions with functions 



with hermite polynomials Ha(x) and sum over all a (in d dimensions one 
has to introduce products of d such functions 

< =%— % (17) 
and to sum over all of them) then one obtains for the term 
Hkl,k2,...k[,k'2,...'^l:J^'^k2^^^'^><''l'^k'2^^^ eHmiuation factor 

rk^,k2,-K,k'2,... = ~ (1^) 

i,j i,j ij 

r{k,k') = exp(-^^(A:-A:')Vfco)- (19) 

Use is made of the generating function (19.12.14) in 0. Here p and kQ are 
appropriately chosen parameters. This function r has the following nice 
properties: If the momentum A; of a creation operator and the momentum 
k' of an annihilation operator are equal, then the function r is the same if 
we remove these two operators. Therefore multiplication of a term with the 
number operator c|.Cfc will not change r. Consequently if all the momenta k 
of and k' of c are pairwise equal, then the factor r vanishes. Indeed these 
matrix-elements to the hamiltonian are considered to be diagonal. If the 
differences between all momenta are large in comparison to ko, then only 
the contributions for i = j are left, so that r approaches the sum of the 
number of creation and annihilation operators. If the number of creation 
operators is different from the number of annihilation operators, then for 
nearly equal momenta the factor r becomes the square of the difference of 
the number of creation and annihilation operators. 

If a system is to undergo a superconducting transition then one should 
keep pairs of creation operators with momenta k and —k, similarly for anni- 
hilation operators. In this case it is more appropriate to sum only over the 
contributions with odd a in ([T^). In d dimensions one keeps those with odd 
ai + 02 + ...Oil- For these odd a the contribution + U-fc vanishes, so that 
multiplication by c|.cl^ or CkC-k does not alter the factor r. If a matrix- 
element consists only of such factors and factors c^Cfc, then its r vanishes. 
With this choice one obtains the elimination factor (|l^) with 

r{k, k') = \ (^exp(-^(A; - l^flkl) - exp(-^(A; + k!f jk^^ . (20) 
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It should be mentioned, that G need not necessarily be a quadratic form 
in the Hamiltonian. For a general differentiable G one obtains 



r]=[H,R] (23) 



Choosing now 

we obtain 
with 

The choice 

G = -J2kHkk (25) 
orders the states with increasing energieJil. One obtains 

Rji = -jSji, Vji = - (26) 

Other choices of G linear in H are the introd iictior i of tjji = sign(i — j)Hji 
which also orders the states in increasingjordeiil3E2l or choosing R being the 
double-occupancy in the Hubbard modeled. 

Note that the introduction of G gives a straighforward way to obtain 
an rj which minimizes G. However it is not the only way to do—this. A 
counterexample is the modification used for the spin-boson modeled. 



3. SECOND-ORDER PERTURBATON THEORY 

We now use this scheme in order to determine the flow equations 

^ = m,Hil)] = mVrlH] = [[T + V,Vr],T + V] (27) 
up to second-order in the interaction V 

H = T + V (28) 
T = Y.^, : c\,c,s : (29) 

q,s 

^ = ^ yiki,k2,qi,q2) ■■cl^siCq,sicl^s2^q^s2 ■■ ■ (30) 

ki,qi,si,k2,q2,S2 
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where Q is the volume of the system. In first order 

dFi(0 



dl 

one obtains 



[[T,Vir],T] (31) 



= -Ae[ki,k2,qi,q2) r{ki,k2,qi,q2)Vi[ki,k2,qi,q2;l) 

(32) 

with 

Ae(A;i, A;2, gi, 92) = + e^^ - e^^ - Sq^ (33) 

and 

Vi{ki,k2,qi,q2;l) = e-^^('=i'*^^'«i'^^)'^('=i''=^'^i'«^)V(A:i,fc2,gi,g2). (34) 

Thus the decay of the off-diagonal terms does not only contain the factor 
(Ae)^ in the exponential, but also the elimination factor r. The one-particle 
contribution T does not change in first order in V. We will neglect the 
second-order contribution to T. 
In second order in V we have 

- [[T,V2r],T] = [[T,V,,],Vi] + [[Vi,V,,],T] (35) 

Prom this equation we obtain 

(^ + Ae{ki,k2,qi,q2fr{ki,k2,qi,q2))V2iki,k2,qi,q2;l) = 
^ J2 /(^i'^2,Pi,P2;Pi,P2, 51,^2; 0(1 - - np^) 

P1,P2 

2 ^ 

+ f{ki,Pi,qi,P2;k2,P2,q2,Pi;l){np^ -"-^2) 

~7l Yl fi^i^Pi^P2,<li;k2,P2,q2,Pi;l)inpi - rip^) 

Pl,P2 

~7j fiki^Pi^<li^P2;k2,P2,Pi,q2;l)inpj^-np2) 

Pl,P2 

m /(^i'^'i'P2,g2;A;2,P2,Pi,gi;0('^j'i - ^P2) (36) 



n 

P1,P2 



with 



f{{a}-{b};l)=Vi{{a}-l)Vi{{b};l) 
X {r{{a}){2Ae{{a}) + Ae{{b})) - ri{b}){Ae{{a}) + 2Ae{{b}))) . (37) 
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For those momenta which obey r{ki,k2,qi,q2) = we obtain 

V2{ki,k2,qi,q2;oo) = 
- F{ki,k2,pi,P2;pi,P2,qi,q2){l - rip^ - Up^) 

Pl,P2 

+7T 2^ Fiki,Pi,qi,P2;k2,P2,q2,Pi){np 

Pl,P2 



H F'{ki,pi,p2,qi;k2,P2,q2,Pi)inp^ - rip 

Pl,P2 

F{ki,pi,qi,p2;k2,P2,Pi,q2){np^ - np 



Pl,P2 



7^ ^iki, Pi, P2,q2;k2,P2, Pi, qi)inpi - ^pi) (38) 



Pl,P2 



with 



Fi{a};{b})=Vi{a})V{{b}) 
, r{{a}){2Ae{{a}) + Ae{{b})) - r({fe})(A6({a}) + 2Ae{{b})) 
r({a})Ae({a})2 + r({6})A6({6})2 



(39) 



Since the initial condition V{ki,k2,qi,q2) = U6k-^^k2,qi+q2 °f Hubbard 
model as well as r and Ae are invariant upon exchange of the last two 
arguments and independently against exchange of the first two arguments, 
the same holds for Vi. Consequently F obeys the symmetry relations 

F(ai, 02,03,04; {6}) = F(a2,oi,03,04;{6}) = F(ai, 02, 04, 03; {6}) (40) 
F{{a};bi,b2,b3,b4) = F{{a};b2,bi,b3,b4) = F{{a};bi,b2,b4,b3). (41) 

Then the second, third and fourth term in the expression (R^) for V2 cancel. 



3.1. Effective Interaction in various Channels 

We wish now to retain all terms in the interaction which seem to be 
important for phase-transitions in the Hubbard model. In particular we are 
interested to keep those relevant for antiferromagnetism and superconduc- 
tivity that is terms of the form cjj_,_g^CfccJcg+g(, and cl.c^_i^C-qCq. We use the 
elimination factor r in the form 

r{ki,k2,qi,q2) = E (^"(^1) + ^"(^2) - v"{qi) - v^{q2)f . (42) 



stability Analysis of Hubbard Model 



In order to keep the interaction terms relevant for antiferromagnetism we 
choose 

^;"(g) = ^;-(g + go). (43) 

Further, in order to keep the pair-interaction we require 

= -v'^i-q). (44) 

With this choice the contributions hsted in the first column of table 1 will 
be kept. Obviously, Vq contains the effective interaction between pairs of 
electrons of total momentum 0, Va contains the effective antiferromagnetic 
interaction with staggered wave- vector qq. Vc and Va contain also charge 
waves. The contributions Vh and Vp contain the contributions kept in the 
Fermi-liquid picture. With the choice (^) and ( |4^ ) also the terms Vy sur- 
vive, which describe the pair-interaction of electron pairs with total momen- 
tum go- 

With this choice of the elimination function r it turns out, that in all 
cases in the table we obtain r({a}) = r{{b}). Then these factors cancel 
altogether and the expression for F({a}; {b}) reduces to 



V 


W 


{a} 


{b} 


V2B{k,q) = V2(k-k,q,-q) 
V2B = Wb 




Wb 


(fc,-fc,pi,P2) 

(fc,pi,p2,-g) 


ipi,P2,q,-q) 
(-k,p2,pi,q) 


V2H(k,q) = V2(k,q,k,q) 
V2H=Wh+Wf 


Wp 
Wh 


(k,q,p\,p2) 
(k,pi,P2,q) 


{pi,P2,k,q) 
iq,P2,Pl,k) 


V2F{k,q)=V2{k,q,q,k) 
V2F = Wp 


Wp 



{k,q,Pi,P2) 
{k,pi,p2,k) 


{pi,P2,q,k) 

{q,P2,pi,q) 


V2A{k,q)=V2ik,q+qo,q,k+qo) 
V2A = Wa 


Wa 



{k,q+qo,pi,P2) 
{k,Pi,P2,k+qo) 


{pi,P2,q,k+qo) 

{q+qo,P2,pi,q) 


V2C (k,q) = V2 {k,q+qo,k+qo,q) 
V2C=Wc + Wa 


Wa 
Wc 


{k,q+qo,pi,p2) 
{k,pi,p2,q) 


{pi,P2,k+qo,q) 
{q+qo,P2 ,pi,k+qo) 


V2Y{k,q) = V2ik,qo-k,q,qo-q) 
V2Y = Wy 



Wy 


(k,qo-k,pi,p2) 
(.k,pi,P2,qo-q) 


ipi,P2,q,qo-q) 

{qo-k,P2,Pi,q) 



Table 1. Matrix-elements in second order perturbation theory. In the first 
line the data of the first term in (^8|) are given, in the second line those of 
the fifth term. In a number of cases r = 0. These terms do not contribute. 
They are indicated by 11^ = in the table. 
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This has the big advantage that we need not choose the functions v 
expUcitely, which prevents us from any arbitrariness in the evaluation. 
The terms in the upper hues can be rewritten 

^F.A = F{ki,k2,pi,P2;pi,P2,qi,q2){l-np 



P1.P2 



P2 



P2 



_|_ (■ £"-^'' ^2 



Pl+P2,S 



(46) 



usmg 



and 



Ae{{a}) 

Mm 



^Pi 



^P2 



I 

~r ^pi ~ e 



(47) 
(48) 







S 


Wa 


Cfc + fig Cfe + eg 
£fc + Sg+go ~^ ^k+qo 


k + q 
k + q + qo 



Table 2. e", e^ and s in eq. (M) 



The terms in the lower line can be written 
W^B,H,c,Y = -7^^ F{ki,pi,p2,q2\k2,P2,Pi,qi){np^-np^) 



with 



U 



P1.P2 
2 



E 



(np-^ ) (e 



Pi + J ) 



(epi + 2 



2~ 



(49) 



and 



Ae({a}) 
Ae({6}) 



e + Epj^ ep2 



Cp2 "^pi 









s 




Cfc — eg 


eg — efc 


A; + g 




Cfc — eg 


efc — eg 


k — q 




Cfc — eg 


efc+go ~ ^q+qo 


k — q 




efc — eg(,_q, 


^g ~ ^go-fc 


k + q-qo 



(50) 
(51) 



Table 3. e", e^ and s in eq. (|9D. 
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Using e_p = Cp both terms can be written in the form 

^ (e,, + - + i^? ^ ^ 

For the first term Wf,a we have 

7lp2 1T'P21 ^Pl ~ fpi ) (^3) 

whereas for the last term VFb,h,c,y one puts 

-en- (54) 



3.2. Isotropic Model 

Before going for the Hubbard Model let us consider an isotropic electron 
gas without any lattice with one particle energy and interaction 

efc = V{ki,k2,qi,q2) = U6k^+k2,qi+q2- (55) 

The numerical calculations were performed by Kording and are taken partly 
fromil. They yield for the three-dimensional system at zero-temperature 
with dimensionless quantities kp = 10, m = 1, the following p-, d-, and 
/- wave-couplings Wb at the Fermi-surface — 0.1174J7^, — 0.01808J7^, and 
— 0.006241C/^, resp. Thus the couplings decay rapidly with increasing angu- 
lar momentum, but they are attractive. 

In two dimensions the numerical calculation yields as a function of tem- 
perature 



1/T 


p-wave 


d-wave 




coefficient of 




cos{e) 


C/2 cos(20) 


1 


-2.18 • 10"^ 


-2.02 • 10-^ 


2 


-1.10 • 10"^ 


-1.06- IQ-^ 


3 


-7.32 • 10"^ 


-7.14- 10"^ 


4 


-5.50 • 10"^ 


-5.39 • 10"^ 



Table 4. Effective pairing interaction in the two-dimensional isotropic 

model. 

One finds, that the contributions decay approximately proportional to 
T. The same tendency can be seen for the ti^and g-wave. This is in agree- 
ment with the prediction by Feldman et alBa, that at zero-temperature in 
second order of the coupling only the s-wave part of the pair-interaction is 
changed, but no contribution to higher angular momenta is obtained. 
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4. HUBBARD MODEL 

Explicit calculations were performed for the Hubbard model on a square 
lattice described by the Hamiltonian 

H = -t J2 4's(^rs + U^{nr^ -^){nri-^), (56) 

{r',r),s i~ 

where the hopping-terms are summed over nearest neighbors. 



4.1. Symmetries 

As long as ego_fc = — we find a number of identities. Let us begin 
with ^ 

If one replaces under the sum pi by pi + and p2 by p2 — Qo and considers 
only the 1 in the first paranthesis, one observes that this contribution changes 
into its negative. Therefore we may skip the 1. Moreover the expression is 
symmetric with respect to the exchange of pi and p2 ■ Thus we may write 

If we now replace p2 by P2 + Qo, then we obtain 

Similarly the expressions for Wb and VFc can be rewritten 

Pl,P2 Pi P2'^ 1^ 

- i^P.-^P.f^i^^-^y''-"'^'-'- ^ 
Prom these expressions it is apparent, that 

WBiKq) = Wc{K-q) = -WA{k,q) (62) 

A second similar relation is obtained for Wh and Wy- In these cases 
e» = and 

WYi{Kq) = WY{k,qo-q) (63) 
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hold. If moreover = (half- filling) , then by replacing pi to qo — pi one 
obtains 

WF{k,q) = -WY{k,q). (64) 
Another class of symmetries is the following 

WA,c,Y{k, q) = Wa,c,y(/c + go, g + go)- (65) 

In the first two cases e"" and interchange. Thus the equation holds. In the 
last case and —e^ have to be exchanged. If we simultaneously exchange 
the momenta pi and —p2 in (^9|), then also the last equality is obtained. 



4.2. Free Energy in Molecular Field Approximation 

After having calculated the effective interaction we are left with a prob- 
lem which is much simpler then the original one. The reason is that the 
original two-particle interaction is a function of three independent momenta. 
Here, however, the interaction depends only on two independent momenta. 
Thus the fluctuations, which make calculations for the original interaction 
difficult are basically eliminated. It is now sufficient to find the minimum of 
the free energy within a Hartree-Fock-Bogoliubov approximation. In prac- 
tice we will start from the symmetric state and investigate, whether this 
state is stable against fluctuations of the order-parameters. Thus we have 
to expand the free energy 

F = {H)- TS, (66) 

in the fluctuations around the symmetric state (S* denotes the entropy). For 
this purpose we represent {H) as a function of the expectations 

{ala]) = Am, Am = -Ai^, (aiak) = A*m (67) 
(alai) = nl6M + m, m = J^ik- (68) 

We assume A and to be small quantities and expand the entropy up to 
second order in these quantities. After some calculation one obtains 

4')), (69) 



k 

- ^-jT.\^kk'?f{(i{el-l^),fi{el.- 

kk' 

- ^^|A,fc,|V(/3(e^-/x),/3(/^- 

kk' 
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with 

S' = -fcB^(nOlnnO + (l-nO)ln(l-nO)), (70) 
k 

nx,y) = ^:^ =^^(e- + l)(e^ + l). (71) 

In total we have to consider four channels, two particle-hole and two particle- 
particle-channels, since the total momentum can be and qq. In all cases 
we have to distinguish between singlet- and triplet-excitations. 



4-2.1. Particle-particle channel, q = 0, interaction Vb 
In this channel we consider the expectation values 

{4,s/-kJ = e,„s,Af + {ea-U,,,Ar* (72) 
(c_fc,,,Cfe,,i) = es,,sAk + i(^'"^)s2,sAk" (73) 

with € = a^. They obey the symmetries 

A| = Ai„ Al- = -A%. (74) 

Due to these symmetries s- and d-wave pairing is of singlet type, p-wave 
pairing of triplet type. The interaction reads 

^B = ^ ^ VBik,q)cl^Cgscl^^,c-gs', (75) 

k,q,s,s' 

which yields the energy 

= ^ E ^B(/c,g)(e,,,Ar + (ea"),,,Ar)(e,,,A^ + (a/^e),,,Af) 

k,q,s,s' 

= ^T.Mk,q){A%*Al + ArAl"). (76) 

k,q 



and the entropy 



Sb = -kBj2fB{k){Al*Al + Al-*Al-), (77) 

k 

fB{k) = f{fi{ek-ii),l3{ii-ek)). (78) 
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4-2.2. Particle-hole channel, q = 0, interactions Vr and Vp 
In this channel the expectation values 

{clcks') = Ss,Anl + 4) + (79) 
are considered. They obey the symmetries 

4* = 4, yr = ^T (80) 

The corresponding part of the interaction reads 

^HF = ^ II ^H(A:,g) : 4,Cfc,cJ^,Cg,/ : 

+ ^ II VF(A;,g) : 4^Cg,cJ^,Cfcy : . (81) 

k,q,s,s' 

Its contribution to the energy is 
2 v-^ 

^HF = ^2^VH(fe,g)z^feZ^^ 

- ^ E ^F(fc,g)(<^.,3'^I + <3'^^r)('^a',.^^ + ^f..^^f) 

= ^E(2^H(fc,9)-^F(fc,g))^I^^-^E^F(fc,g)4"^r' (82) 

fe,q fe,g 



and its entropy reads 



5hf = -fcBE/H(^)((^D' + (4")') (83) 
/H(fc) = /(/3(eik-/x),/?(efe-/x)). (84) 



4-2.3. Particle-particle channel, momentum qo, interaction Vy 

In this channel the expectation values 

= ess'A%* + {ea^),,Ar (85) 
{cgo-ks'Cks) = e,/,A| + ((7"e)y,A*fc°^ (86) 

are considered with the symmetries 

Al^_, = Al A*^_, = -Ar. (87) 
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The interaction reads 

k,q,s,s' 

and its contribution to the energy 

k,q,s,s' 

= ^E^Y(A^,<?)(ArA^ + ArAr). (89) 

k,q 

Let us now divide the Brillouin zone into two halves. The contributions 
to the sum with \k\ < l^o ~ ^| are denoted by a prime '. Then the energy 
reads by use of eq. (|6^) 

k,q 

+^ E go - 'z)(Ar A^ - Ar A^") 

k,q 



^(EW(fc,g) + W(fc,go-g))ArA^ 

k,q 



+^E - VYik,qo - g))Ar A*", (90) 

k,q 

and the corresponding entropy 

= "Y E /Y(^)(e-'Ar + (6a"),,,Ar)(es'sA| + (^^6),,,Af) 

k,s,s' 

= -kBY.fY{k){AfAl + A'rAi») 

k 

= -2A;BE/Y(fc)(ArA| + ArAr), (91) 

k 

fY{k) = f{P{ek-^^),P{^l-eg,+k). (92) 



4-2.4- Particle-hole channel, momentum qo, interactions Va and Vc 
Here the expectation values 

{clck+q,s') = Sss'iyl + aZ,ul'' (93) 
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are considered, which obey the symmetries 

^U,,=^k ^l%o = ^'r- (94) 
The interaction reads 



k,q,s,s' 

2Q 



+ ^ 5Z ycik,q)cl^Ck+g^scl+g,s'(^qs'- (95) 



k,q,s,s' 

Its contribution to the energy is 

^AC = E ^A(A;,g)(5,,,z.| + <,z.r)(5,,,z.- + a,%z.f*) 

k,q,s,s' 

+ ^ E ^C(^,'Z)^K (96) 

k,q,s,s* 

= ^ E (2^c(fe, g) - V^a(A:, q)HiyT - ^ E ^a(A;, g)^"^^*- 

fc,(jr k,q 



Summation over half of the Brillouin zone yields with eq. (|6 

= ^E(2^c(fc,'?)-i^A(fc,g))(^^r+^r^P 

k,q 

+ ^ E (2^c(fc, <? + go) - i^A(fc, q + go))(^i^? + z^ro 

k,q 

- ^El^A(A:,g)(4"^r + ^fc"*^f) 

- ^E^A(fe,g + 9o)(4"< + 4"*^r)- (97) 
We now separate the real and the imaginary part of i' 

T^k = ;^(^'^ + ^^+90 ^ ^^'^^'^ ~ ^^'^^ *^98) 

i^fcl^g + i^fci^g = T^kVq - Vk^q (99) 
Vk^l + vlVq = Vkl^q + VkVq. (100) 

Then the energy reads 



^AC = ^E(2^c(A;,g)-^A(A:,g) + 2yc(A:,g + 'Zo)-^A(A:,(? + go))i^fci^g 

k,q 
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+ E (2^c(fe, q) - VAik, q) - 2Vc{k, q + qo) + Va(A;, q + qo)M^l 
k,q 

k,q 

- ^E(^A(fc>5)-^A(A:,g + go))5r5r- (101) 

k,q 

The corresponding entropy is 

k,s,s' 
k 

= -kB E fc{k)mf + ipif + ipr? + {^\"f) 

k 

fcik) = /(/3(efe-/i),/3(efe+,o-/^))- (102) 

4.3. Variation of the Free Energy 

The expression for the free energy which is bihnear in i/ and A, resp., 
has to be checked with respect to its stabihty. That is, as soon as some v or 
A different from zero yields a lower free energy than for the symmetric state 
for which all u and A vanish, then the symmetric state is unstable and the 
system will approach a symmetry broken state. This is the indication for a 
phase transition. The expression for the free energy has the form 

= ^J2((^U)il + jVk,g)AlA, + J2f,AlAk (103) 

k,q k 
E {j^k^l + ^jfBk^l + ^k,^ V^^fe \[Tq^q (104) 



k,q 



with 



(5t ^ _ pm,q 



Ako = — Bka = ^-^=^- (105) 

A similar bilinear contribution is obtained for u instead of A, which is han- 
dled in the same way. Here the factor U"^ /t is extracted from the matrix 
elements V^^q- These matrix elements and the entropy coefficients depend 
on (5t and (iji. The same is true for the coefficients A}^^q and B^^q. 
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The calculation is performed for the various representations under the 
group C4 = 4mm. The representations of the even-parity states are one-di- 
mensional. We denote them hy s+ = si, s- = g = s^y(^^2_y2-^, d+ = d^2_y2, 
d- = dxy The odd-parity representation is two-dimensional, here simply 
denoted by p. Moreover in the channels VX, Vc, and Vy we can distinguish 
between eigen solutions Vk = ^^'k+qo-, and = ztAfc+gg, resp. In total the 
calculation is performed for 45 channels. 

For each channel the eigenvalues A of ^A+{^-)'^B have to be determined. 
Whenever the lowest eigenvalue (i.e. the most negative) equals —1, then a 
critical {U/t)c is reached. For the s+ representation (in case of Va, Vc, 
Vy only for = +Vk+qo, = +^k+qo) one has to find the solution by 
iterating the eigenvalue equation as a function of U /t. In the case of all 
other representations the j4-term does not contribute. Therefore then we 
determine the lowest eigenvalue A of B, and obtain {U/t)c = A. 



5. NUMERICAL RESULTS 
5.1. Numerics 

In a first step we determine the matrix elements as given above on a 
grid of 2nh x 2nh lattice points in the Brillouin zone. The resultsE3 presented 
below are calculated for = 16. The calculation of the matrix-elements 
takes most of the computer time. It increases with the sixth power of nt. 

It is obvious that in a number of cases the denominator becomes ex- 
tremely small or even vanishes. We remember, that the sums (^) are of the 
form 

E;fe- (106) 

p ^p^ ^ 

where e is independent of p. In particular if e vanishes then we have basically 
a main value integral 

/ d^p^ (107) 

where Zp can vanish. In order to avoid, that we sum up terms close to a 
zero of Zp, which would yield eratic contributions we proceed as follows: We 
average over several points in momentum-space in the vicinity of p. More 
precisely if we have to sum 

^'"^'^^'^^ ^k+q,pi+p2, (108) 
pi,P2 ''^k,q,pi,p2 

then we calculate z and n not only for the given {k, q,pi,p2), but also for 
{k + S,q,pi+6,p2), {k + 6,q,pi,p2 + 6), {k,q + 6,pi + 6,p2), {k,q+6,pi,p2 + 6), 
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{k + 6,q — 6,pi,p2), {k,q,pi + 5,p2 — 6) for four different 6s. They are 
6 = (±^,0) and 6 = (0, it^), where 6p is the momentum-spacing for the 
evaluation of the sums. Then we attribute to these the weights gn and to 
{k, q,pi,p2) the corresponding weight n, which means that we use 

pi,P2 ^k,q,pi,p2 + 9 2^{d} n..±s.. 

with some choice of g. Since n is never negative we obtain nearly always a 
positive denominator. An exception are the Vn and l^-terms which obey 
k = q = p and 2k = qo modulo reciprocal lattice vector, and the V^-terms 
with p = q — qo = —k and 2k = qo (modulo reciprocal lattice vector). In 
these cases all ns vanish and we replace the term by the average of the four 
terms at p ± {5p, 0) and p ± (0, 6p). 

Since the matrix-elements are calculated by summing over a grid in the 
Brillouin zone, the question arises, how reliable are the results? Obviously 
at low temperatures the occupation numbers n vary rapidly. So the result 
may depend on the points in our grid. There are basically two criteria: 
(i) In the channels V2F, 1^2H, and V2Y the denominator in (106) contains an e 



which vanishes identically. Despite of the procedure described above one has 
to assume, that the calculated matrix-elements are less precise than those 
for Vb, Va, and Vq. 

(ii) The increase of the entropy factors / determines the extension in the 
wave- vector space contributing essentially to the order-parameter. Denoting 
X = P{e — fj.) one observes, that /b = f{x,—x) increases like 2\x\, whereas 
/h = f{x,x) increases like e'^' for large Thus the contributing phase 
space for Vq is larger than for Vh,f- The other interactions lie in between. 
Indeed it turns out, that the results for the particle-hole channels coming 
from Va^F show deviations between different nhS at low temperatures, that 
is below O.OSt. 



5.2. Results 

In the following figures we show the critical values of the lowest lying 
(U /t)c as function of T/t for different values of /x/t. We have chosen g = 0.01 
in (|109| ). Calculations for 51 = 1/24 instead yield practically no difference. 

In the following table we summarize the frequently observed low-lying 
instabilities. In the second column we charactarize the symmetries: ph 
stands for particle-hole, pp for particle-particle, si for singlet, tr for triplet, 
q± for i^k+qo = ^^k, and ^k+qo = iA^, otherwise 0. If there are instabili- 
ties in another channel with an (U/t)c less than that for the Pomeranchuk- 




Fig. 1. {U/t)c at n/t = 0. Doping x = 0. Also shown: At T/t = 0.01 
symmetry {ph tr 5+) (ferromagnetism, disappears for larger lattices). 

instability with (s_ = g) and less than 10, then the lowest one of those is 
also shown. 





channel 


color 


Antiferromagnetism 
Pomeranchuk instability 
Band splitting 
Superconductivity 
Flux phases 
Pomeranchuk instability 


ph tr q+ s+ 
ph si d+ 
ph si q- p 
pp si d+ 
ph si/tr g_ d+ 
ph si S- = g 





Table 5. Channels with low-lying (U/t)c. 



As expected at half- filling (/i/t = 0) the lowest lying is the antifer- 
romagnetic one. It turns out, that the second order contribution to the 
antiferromagnetic channel suppresses antiferromagnetism (at least in the s- 
channel). Therefore antiferromagnetism disappears at larger values of U/t. 
This is expected since in the strong-coupling limit the magnetic interaction 
vanishes with J = At'^/U in the t-J-model. As we leave half-filling the 
antiferromagnetism becomes weaker and disappears at some doping. This 
observation is remarkable, since although we work with a weak-coupling 
calculation, we obtain results expected at strong couplings. 

The next instability is a Pomeranchuk-instability. It has recently been 
observed from RG-calculations by Halboth and Metzneill3. It should be 
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Fig. 3. {U/t)c at ^/t = 1/8. Doping varies from x = 0.0807 to x = 0.0519. 
Also shown: At T/t = 0.01 symmetry {ph tr p). From T/t = 0.02 to 0.04 
with symmetry (pp si S- = g). 
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10- 
8- 
6- 
4- 
2- 



0.00 0.05 0.10 0.15 0.20 0.25 0.30 

Tit 

Fig. 4. {U/t)c at ii/t = 1/4. Doping varies from x = 0.1159 to x = 0.1030. 

mentioned, that our calculations become less reliable for the p/iO-sector at 
very low temperatures. Thus we do not know, whether (except from the 
antiferromagnetic instability) it is the lowest lying one as T approaches zero. 
Depending on the choice of nh the slope of the curve remains positive at low 
temperature or it becomes negative. 

The next instability is a particle-hole instability of singlet type with g_ 
and p-type symmetry. It corresponds to an strengthening and weakening of 
the hopping matrix-elements at alterning bonds along the x- or y-direction 
or both superimposed. This instability leads to a splitting into two bands. 
Within the present approximation they are not separated by a gap. One 
may speculate, that in higher orders they develop into the two Hubbard 
bands. We observe, however, that at sufficiently low temperatures either 
antiferromagnetism or d-wave superconductivity yield a lower critical (U/t). 

Then the superconducting dx2_y2 instability follows, which is much 
stronger than superconducting instabilities of other symmetries. 



symmetry 


T 


= O.lt 


T 


= 0.03t 


s+ 


12.51 .. 


.-0.10 


25.90 . 


..-0.39 


s-= g 


0.00 . . 


.-0.71 


0.00 . 


..-1.75 


d+ 


0.01 .. 


.-2.49 


0.01 . 


..-5.49 


d_ 


0.69 . . 


.-0.04 


2.36 . 


..-0.06 


P 


1.35 .. 


.-0.32 


2.74 . 


..-0.97 



Table 6. The range of the eigenvalues A at half-fillingEj multiplied by 100. 
Here only the second-order contribution is taken into account for s+. 
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Finally we observe a flux-phase instability, which has been discussed by 
KotUai0, and by Affleck and Marstonlll. It has been recently discussed by 
Chakravarty et al. as d-density wave-ordert3. The singlet and triplet eigen- 
values are degenerate for even parities (apart from [q^ s+)), since according 
to table 1 and eqs. (|6^Jl0l| ) the singlet interaction 2V2C — ^2A = 2Wc+Wa = 
— Wa equals the triplet interaction — V2A- At /i = it is moreover degenerate 
with the superconducting instability, since V2B = Wb = —Wa = —V2A and 
for vanishing fi also the entropy factors /b = /c agree. 

At even higher values another Pomeranchuk-instability appears with 
S- = g wave-character. A few other low-lying instabilities observed only 
at special parameters are shown in the figures in black. The ferromagnetic 
instability at T = O.Oli and = may be a remnant of the Nagaoka- 
ferromagnetismii, although it is expected at larger values of U /t. For = 
24 it does no longer appear in our calculation even at T/t = 0.01. 

One phase may suppress another phase. To which extend two order- 
parameters can coexist with each other is another question, which should be 
investigated in the future. 



5.3. Conclusion 

In these calculations we obtain many critical couplings (U/t)c of the 
order of 5. Since this is not a small number, calculations which do not 
handle the coupling perturbatively, should be performed. Nevertheless it 
becomes clear that this type of calculation gives a good estimate of the 
most important instabilities. We have found all commonly dicussed types of 
order. Moreover, a p-wave instability which yields band-splitting, appears. 
Although our calculation is performed for weak coupling, a number of effects 
normally obtained in the strong-coupling limit are reproduced reasonably 
well by means of the flow equations. 
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